home *** CD-ROM | disk | FTP | other *** search
/ SGI Hot Mix 17 / Hot Mix 17.iso / HM17_SGI / research / lib / crvlength.pro < prev    next >
Text File  |  1997-07-08  |  3KB  |  86 lines
  1. ; $Id: crvlength.pro,v 1.4 1997/01/15 03:11:50 ali Exp $
  2. ;
  3. ; Copyright (c) 1996-1997, Research Systems, Inc.  All rights reserved.
  4. ;       Unauthorized reproduction prohibited.
  5. ;+
  6. ; NAME:
  7. ;       CRVLENGTH
  8. ;
  9. ; PURPOSE:
  10. ;       This function computes the length of a curve with a tabular
  11. ;       representation, y(i) = F(x(i)). 
  12. ;
  13. ; CATEGORY:
  14. ;       Numerical Analysis
  15. ;
  16. ; CALLING SEQUENCE:
  17. ;       Result = Crvlength(X, Y)
  18. ;
  19. ; INPUTS:
  20. ;       X:    An N-element vector (N >= 3) of type float or double. These 
  21. ;             values must be specified in ascending order. Duplicate x values 
  22. ;             will result in a warning message.
  23. ;
  24. ;       Y:    An N-element vector of type float or double.
  25. ;
  26. ; KEYWORD PARAMETERS:
  27. ;       DOUBLE:  If set to a non-zero value, computations are done in
  28. ;                double precision arithmetic.
  29. ;
  30. ; RESTRICTIONS:
  31. ;       Data that is highly oscillatory requires a sufficient number
  32. ;       of samples for an accurate curve length computation.
  33. ;
  34. ; EXAMPLE:
  35. ;       Define a 21-element vector of X-values.
  36. ;         x = [-2.00, -1.50, -1.00, -0.50, 0.00, 0.50, 1.00, 1.50, 2.00, $
  37. ;               2.50,  3.00,  3.50,  4.00, 4.50, 5.00, 5.50, 6.00, 6.50, $
  38. ;               7.00,  7.50,  8.00]
  39. ;       Define a 21-element vector of Y-values.
  40. ;         y = [-2.99, -2.37, -1.64, -0.84, 0.00, 0.84, 1.64, 2.37, 2.99, $
  41. ;               3.48,  3.86,  4.14,  4.33, 4.49, 4.65, 4.85, 5.13, 5.51, $
  42. ;               6.02,  6.64,  7.37]
  43. ;       Compute the length of the curve.
  44. ;         result = CRVLENGTH(x, y)
  45. ;       The result should be:
  46. ;         14.8115
  47. ;
  48. ; MODIFICATION HISTORY:
  49. ;       Written by:  GGS, RSI, March 1996
  50. ;-
  51.  
  52. FUNCTION CrvLength, X, Y, Double = Double
  53.  
  54.   ON_ERROR, 2
  55.  
  56.   TypeX = SIZE(X)
  57.   TypeY = SIZE(Y)
  58.  
  59.   ;Check Y data type.
  60.   if TypeY[TypeY[0]+1] ne 4 and TypeY[TypeY[0]+1] ne 5 then $
  61.     MESSAGE, "Y values must be float or double."
  62.  
  63.   ;Check length.
  64.   if TypeX[TypeX[0]+2] lt 3 then $
  65.     MESSAGE, "X and Y arrays must contain 3 or more elements."
  66.   if TypeX[TypeX[0]+2] ne TypeY[TypeY[0]+2] then $
  67.     MESSAGE, "X and Y arrays must have the same number of elements."
  68.  
  69.   ;Check duplicate values.
  70.   if TypeX[TypeX[0]+2] ne N_ELEMENTS(UNIQ(X)) then $
  71.     MESSAGE, "X array contains duplicate points."
  72.  
  73.   ;If the DOUBLE keyword is not set then the internal precision and
  74.   ;result are identical to the type of input.
  75.   if N_ELEMENTS(Double) eq 0 then $
  76.     Double = (TypeX[TypeX[0]+1] eq 5 or TypeY[TypeY[0]+1] eq 5) 
  77.  
  78.   nX = TypeX[TypeX[0]+2]
  79.   Yprime = (SHIFT(Y,-1) - SHIFT(Y,1)) / (SHIFT(x,-1) - SHIFT(x,1) + 0.0)
  80.   Yprime[0] = (-3.0*Y[0] + 4.0*Y[1] - Y[2]) / (x[2] - x[0])
  81.   Yprime[nX-1] = (3.0*Y[nX-1] - 4.0*Y[nX-2] + Y[nX-3]) / (x[nX-1] - x[nX-3])
  82.  
  83.   RETURN, INT_TABULATED(X, SQRT(1 + Yprime^2), Double = Double)
  84.  
  85. END
  86.